复杂模型事后解析方法

内在可解释机器学习模型本身虽然具有较好的解释性, 但其在精度上依旧无法与复杂模型相媲美。在某些场景下, 如果既想使用复杂模型, 又想要模型的结果具有解释性, 那么我们可以采用事后解析方法, 弥补复杂模型难以被解释的不足。 生物统计中常用的敏感性分析、数理统计中经典的 PDP (Partial Dependence Plot, 部分依赖图) 和 ALE (Accumulated Local Effect, 累积局部效应）、合作博恋论中的公平分配贡献, 以及近几年流行的新方法, 如代理模型等, 这些都可以用来做事后解析。

不得不使用黑盒模型时, 科学家会在原有模型的基础上, 运用各种事后解析方法, 试图打开复杂模型的黑箱子, 让模型使用者更好地理解箱子内部的运作机理。Marco Tulio Riberio、 Sameer Singh 和 Carlos Guestrin 于 2016 年提出了 LIME(Local Interpretable Model-Agnostic Explanation, 与模型无关的局部解释算法)的概念, 即用代理模型的思想对单个样本的模型预测值进行解释。Scott M. Lunberg 和 Su-In Lee 于 2017 年在 NIPS 会议 (神经信息处理系统大会) 上提出 SHAP (ShapleyAdditive Explanation) 的概念, 引人合作博恋论中的 Shapley Value 来计算样本中每个特征的贡献值, 并以此作为局部解释。

部分依赖图

在机器学习的发展过程中, 挖掘数据之间线性关系的线性回归模型, 在关系更为复杂的数据集上, 其精度明显低于逐渐发展起来的现代模型, 如支持向量机、随机森林、梯度增强算法, 等等。然而, 随着一系列现代黑盒模型在工业、金融、医疗等领域的部署和应用, 这些黑盒模型在提升精度的同时, 也导致了结果难以解释的问题：人们看不到模型内部的运作机 理, 对于模型给出的预测结果, 不能进行有效的阐释。

部分依赖图于 2001 年由 Friedman 提出, 是一种通过可视化图像方式挖掘数据之间关系的工具, 它不仅可以直观地同译預测结果与特征数据之间的关系, 而且可以为模型诊断提供一种反馈机制。

部分依赖函数

假设在银行个人客户申请信用评分问题中, 在完成了从数据清洗到模型构建的全部流程之后，我们得到了一份表现不错的客户末来一段时间内发生违约行为的概率评分结果。这时公司的业务人员除了得到这样一份评分结果之外, 还想要获得一些能显著区分客户末来发生违约概率大小的指标, 以及如何使用这种指标构建一些规则这样就可以在业务层面的筛选环节完成更优质客户的筛选工作，降低时间成本与经济成本。

部分依赖图所勾勒的预测结果与数据之间的关系, 就可以为该业务人员提供他所要的指导信息。具体而言, 部分依赖图所描绘的是, 一或两个目标特征对于模型预测结果的边际影响, 这种边际影响在图像上所反映出来的.就是目标特征对于预测结果的重要程度。

那么, 如何才能反映这种关系呢? 部分依赖函数提供了如下的实现方法指导： $$ f\left(X_S\right)=E_{X_S}\left[f\left(X_S, X_C\right)\right]=\int f\left(X_S, X_C\right) p\left(X_C\right) \mathrm{d} X_C $$ 其中, $f$ 是已经构建好的模型，其输出的是预测结果， $X_S$ 与 $X_C$ 分别代表我们希望探索关系的目标特征与非目标特征。从 式子最右边的一项我们可以发现, 部分依赖函数的本质就是对目标特征 $X_S$ 求边际期望, 若去掉 $p\left(X_C\right)$ 项, 则所求得的就是目标特征 $X_S$ 的边际分布。 对于部分依赖函数更通俗的解释是, 在保证非目标特征与模型预测结果无关的情况下, 找到目标特征与预测结果的关系。

估计方法

部分依赖函数这一数学工具体现的思想非常直观, 但是在实际应用中直接套用这一公式其实是比较难的, 因为很有可能$X_C$ 的大部分取值, 都没有出现在我们有限的数据集中, 从而无法直接计算积分。 在具体实现中, 我们对目标特征的边际期望采用的估计方 式为: $$ f\left(X_S\right) \approx \frac{1}{n} \sum_{i=1}^n f\left(X_S, X_C^{(i)}\right) $$ 其中, $n$ 为数据集中样本的数量, $X_C^{(i)}$ 为非目标特征的第 $i$ 个 样本的取值。

下面就来举例说明。如果有一组数据集中包含人的身高、 年龄、体重三个特征, 预测目标为饭量, 我们想知道体重对于人饭量的影响, 那么体重就是我们的目标特征 $X_S$, 其他的特征就是非目标特征 $X_C$, 根据上式, 我们只需要将数据集中的年龄和体重取平均值 (假设分别记为 $\overline{\mathrm{age}}$ 和 $\overline{\mathrm{height}}$ ), 接下来, 在固定 $X_C=[\overline{\mathrm{age}}, \overline{\text { height }}]$ 的情况下, 依次求得不同身高 取值时模型对饭量的预测结果, 就可以勾勒出身高与饭量的部分依赖图。 这里需要额外提醒的一点是, 在勾勒图像时, 关于身高的取值方式, 我们可以采用数据集中的不同身高取值作为横坐标, 如果数据集中的身高分布偏差较大, 那么我们绘制出的图像可能会与真实关系存在较大的偏差，这里采用的方法一般是，选取数据集中身高的最大值与最小值, 在其中切分出等距的若干个格点, 作为横坐标, 这样可以在很大程度上避免偏差。

通过估计公式，我们可以更进一步发现, 对于部分依赖函数中的 $p\left(X_C\right)$, 假设特征 $X_C$ 共有 $T$ 种取值,若 $X_c^{(t)}$ 这一特定取值在数据集中出现了 $c^{(t)}$ 次, 则 $p\left(X_C^{(t)}\right)=$ $\frac{c^{(c)}}{n}$ 。换言之, 在数据有限的情况下, 部分依赖函数的形式可以变为: $$ \begin{aligned} f\left(X_S\right) & \approx \sum_{t=1}^T f\left(X_S, X_C^{(t)}\right) p\left(X_C^{(t)}\right) \approx \sum_{i=1}^T f\left(X_S, X_C^{(t)}\right) \frac{c^{(t)}}{n} \\ &=\frac{1}{n} \sum_{i=1}^n f\left(X_S, X_C^{(i)}\right) \end{aligned} $$ 所以, 在数据有限的情况下, 若我们用特征值在数据集中 出现的频率作为 $p$ 的估算。

在挖掘了估计方法的内在逻辑之后, 我们还可以更进一步, 我们所掌握的样本数量越多, 对于 $p\left(X_C^{(c)}\right)$ 通过 $\frac{c^{(t)}}{n}$ 的估计就会更加接近真实值, 从而使得部分依赖图表现出的关系更加趋近于真实的关系。

部分依赖图的局限

作为一种可视化的描述方法, 部分依赖图有一个非常明显的局限, 就是其能反映的目标特征数量有限。受制于人类空间思维的局限, 部分依赖图最多只能同时反映两个目标特征与预测结果的关系 (两个目标特征可视化时反映为三维图像), 这样导致的问题是, 当数据集的特征维度非常高时, 我们若想挖掘各特征与预测结果的关系，那么遍历一次将会耗费大量的时间。

而第二个局限, 其实是反映在部分依赖函数中。我们之前提到过, 部分依赖函数本质上是在求解目标特征 $X_S$ 的边际期望, 这其中潜在的假设, 是目标特征与非目标特征之间不具有相关性。若非目标特征与目标特征之间具有相关性，则在函数 中反映出来的关系就变为： $f\left(X_S\right)=E_{X_S}\left[f\left(X_S, X_C\right)\right]=\int f\left(X_S, X_C\left(X_S\right)\right) p\left(X_C\left(X_S\right)\right) \mathrm{d} X_{\mathrm{C}}$ 这样对 $X_C$ 取积分就会无法排除非目标特征对于预测结果的影响, 从而无法求取边际期望。若变量之间存在相关关系, 则会导致在 计算过程中产生过多的无效样本，估计出的值会比实际值偏高。

这里我们继续沿用上节提到的例子来说明这个问题, 即求取人的体重与饭量之间的关系。不过，我们的数据中同时还包含了身高与年龄这两个特征, 很明显, 身高、年茈与体重 是存在巨大相关性的，如果通过式中提到的估算方式去计算体重与饭量之间的部分依赖函数, 那么在等式的右边很可能就会出现 “不可能出现的数据”, 比如, 一个身高 $X_S=$ $130 \mathrm{~cm}$ 的孩子, 在估算过程中会对应到 $X_C=[100 \mathrm{~kg}$40 岁 这样的成年人的数据。

个体条件期望图

前文中提到的部分依赖图实际上是一种对数据的全局解释, 我们确定好非目标特征的值之后，探索单个或两个特征对预测结果的影响，部分依赖图并没有深人探系到单个样本的层面, 因此我们需要借助另一个工具- $\mathrm{ICE}$ 图。

ICE 全称为 Individual Conditional Expectation, 即个体条件期望, 用于捅述个体的特征在变化时其预测值的变化倩况。 个体条件期望图其实就是将关系的变化细化至样本层面的可视化工具。

我们将之前提到的估算方法式稍加改变, 就可以实现对个体条件期望图的勾勒: $$ f\left(X_S^{(j)}\right) \approx \frac{1}{n} \sum_i^n f\left(X_S^{(j)}, X_C^{(i)}\right) \quad j \in\{1,2,3, \cdots, n\} $$ 换言之, 数据集中有多少个样本, 我们就能得到多少个 ICE 曲线, 而对这些曲线取平均值, 得到的就是部分依赖图。

实例演示

本节所采用的例子相较于之前更为具体完善, 该示例将展示如何实现部分依赖图与 ICE 图的可视化。为了使实例清晰易懂, 免于繁杂的数据清洗和处理, 我们在这里使用的是一组样本量较小的数据集, 该数据包含500多个样本与 13 个特征, 13 个特征全部为数值型特征, 且不存在缺失值,房屋价格价格是我们需要预测的 $y$ 。 为了满足在本节一开始提到的 “为黑盒模型提供一定解释性” 的观点, 我们在这里使用 XGBoost 这一黑盒模型对数据进行建模, 将数据集按照 $7: 3$ 分为训练集与测试集, 使用模型默 认参数, 建模完成后, 保存训练好的模型 $f$ 。 由于此数据 集的目标变量 $y$ 为价格, 所以其属于回归问题, 这里选用 XGBRegressor() 进行模型训练, 代码如下:

代码实现

from sklearn.datasets import load_boston
import pandas as pd
from sklearn.model_selection import train_test_split
import xgboost as xgb

data = load_boston()
df  = pd.DataFrame(data['data'],columns=data['feature_names'])
df['target'] = data['target']
X,y = df.iloc[:,:-1],df['target']
X_train,X_test,y_train,y_test = train_test_split(X,y,train_size=0.7)

xgb_reg = xgb.XGBRFRegressor(random_state=1)
xgb_reg.fit(X_train,y_train)
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sn

XC_name = X.drop('NOX',axis=1).columns
XC_mean = X[XC_name].mean().values
XS_min,XS_max = X['NOX'].min(),X['NOX'].max()
XS_grid = np.arange(XS_min,XS_max,(XS_max-XS_min)/1000)
pdp_data = pd.DataFrame(np.tile(XC_mean,(1000,1)))
pdp_data.columns = XC_name
pdp_data.insert(loc=4,column='NOX',value=XS_grid)

下面我们就来绘制部分依赖图, 这里选择的目标特征为NOX, 即我们想要探索这一特征对于价格的影响, 反映到图像上是㤰样的。

回顾之前讲到的估计方法, 首先我们计算除了NOX和价格之外所有特征的均 值, 记为 $\overline{X_C}$, 接下来找到数据中NOX的最大值与最小值, 均匀地划分出格点, 将每一个格点记为 $X_S^{(i)}$, 最后通过训练好的模型, 对所有格点与固定的 $\overline{X_C}$ 所组成的数据进行预测, 即 $f\left(X_S^{(i)}, \overline{X_C}\right)(i=1,2, \cdots)$, 预测结果与格点勾勒出的图像即为部分依赖图。

图中展示的是部分依赖图的结果。在假设区块大小与其他特征都是独立的前提下, 我们可以看到, 在0.6-0.7阶段, 价格随着NOX快速增大, 但随后又减小。

部分依赖图为使用者直观地展示了特征与目标变量的关系, 现在我们已经完成了价格预测结果与NOX这一特征的关系探索, 但是这一关系展示的是NOX的平均效应, 此时若想了解某单一样本的取值情况, 以及NOX与价格的关系, 就需要用到 ICE 图。 绘制 ICE 图时, 我们使用的模型依然是之前训练好的模 型, 这一点无须改变。在部分依赖图中, 我们用目标特征的格 点与非目标特征的均值生成了一组数据集, 但是在绘制 ICE图时, 稍微麻烦的地方就是, 这样的数据集需要生成很多组, 具体的数量是由使用者想要确定的非目标特征取值有多少种来决定的。

在接下来的例子演示中, 假设我们想要探求数据集中的前 500 个样本, 在不同的非目标特征取值(即 $i=1,2 ， \cdots, 500$ 时的 $X_C^{(i)}$ 取值) 情况下, 单个样本对应的 ICE 图, 演示代码如下：

fig = plt.figure()
for i in range(50):
    XC_i = X[XC_name].iloc[i,:].values
    ICE_data = pd.DataFrame(np.tile(XC_i,(1000,1)))
    ICE_data.columns = XC_name
    ICE_data.insert(loc=4,column='NOX',value=XS_grid)
    preds = xgb_reg.predict(ICE_data)
    sn.lineplot(x=XS_grid,y=preds)
plt.xlabel('NOX')
plt.ylabel('Price')
plt.show()

参考资料：